This tutorial demonstrates the BioGSP (Biological Graph Signal Processing) package for analyzing spatial patterns in multiplexed imaging data. We’ll use the built-in toy CODEX dataset to showcase the package’s capabilities including:
# Load required libraries
library(BioGSP)
library(ggplot2)
library(dplyr)
library(tidyr)
library(patchwork)
library(knitr)
# Verify BioGSP package version
cat("BioGSP Package version:", as.character(packageVersion("BioGSP")), "\n")## BioGSP Package version: 1.0.0
BioGSP_Tutorial.R
# Load the built-in toy CODEX dataset
data(codex_toy_data)
# Examine the dataset structure
str(codex_toy_data)## 'data.frame': 18604 obs. of 5 variables:
## $ cellLabel : chr "ROI_0_BCL6BCell_1" "ROI_0_BCL6BCell_2" "ROI_0_BCL6BCell_3" "ROI_0_BCL6BCell_4" ...
## $ Y_cent : num 64.5 71.2 78.3 50.8 52.6 ...
## $ X_cent : num 112 82.6 91.4 92.6 98.8 ...
## $ Annotation5: chr "BCL6- B Cell" "BCL6- B Cell" "BCL6- B Cell" "BCL6- B Cell" ...
## $ ROI_num : chr "ROI_0" "ROI_0" "ROI_0" "ROI_0" ...
## cellLabel Y_cent X_cent Annotation5 ROI_num
## ROI_0.1 ROI_0_BCL6BCell_1 64.52399 111.96648 BCL6- B Cell ROI_0
## ROI_0.2 ROI_0_BCL6BCell_2 71.22437 82.56243 BCL6- B Cell ROI_0
## ROI_0.3 ROI_0_BCL6BCell_3 78.33008 91.44301 BCL6- B Cell ROI_0
## ROI_0.4 ROI_0_BCL6BCell_4 50.83888 92.60675 BCL6- B Cell ROI_0
## ROI_0.5 ROI_0_BCL6BCell_5 52.64693 98.76092 BCL6- B Cell ROI_0
## ROI_0.6 ROI_0_BCL6BCell_6 56.20720 91.39930 BCL6- B Cell ROI_0
# Summary of cell types
cell_summary <- codex_toy_data %>%
group_by(Annotation5) %>%
summarise(
Count = n(),
Mean_X = round(mean(X_cent), 2),
Mean_Y = round(mean(Y_cent), 2),
SD_X = round(sd(X_cent), 2),
SD_Y = round(sd(Y_cent), 2)
) %>%
arrange(desc(Count))
kable(cell_summary, caption = "Summary of Cell Types in Toy CODEX Data")| Annotation5 | Count | Mean_X | Mean_Y | SD_X | SD_Y |
|---|---|---|---|---|---|
| CD4 T | 4092 | 45.64 | 45.52 | 21.49 | 23.61 |
| BCL6- B Cell | 3719 | 55.43 | 45.64 | 22.96 | 20.28 |
| CD8 T | 3346 | 47.38 | 48.53 | 23.91 | 27.34 |
| DC | 2233 | 50.17 | 46.15 | 23.18 | 24.00 |
| CD4 Treg | 1490 | 53.46 | 40.06 | 23.25 | 19.53 |
| M1 | 1490 | 46.59 | 54.36 | 21.19 | 22.62 |
| BCL6+ B Cell | 931 | 36.53 | 51.97 | 21.03 | 23.56 |
| Endothelial | 746 | 51.29 | 45.41 | 22.90 | 19.33 |
| M2 | 370 | 56.67 | 37.24 | 28.21 | 19.59 |
| Myeloid | 186 | 47.60 | 39.69 | 21.09 | 20.63 |
| Other | 1 | 107.12 | 46.12 | NA | NA |
# Focus on BCL6- B Cell and CD4 T for detailed analysis
target_cells <- c("BCL6- B Cell", "CD4 T")
cat("\nTarget cell types for SGWT analysis:\n")##
## Target cell types for SGWT analysis:
## BCL6- B Cell: 3719
## CD4 T cells: 4092
##
## ROI summary:
##
## ROI_0 ROI_1 ROI_10 ROI_11 ROI_12 ROI_13 ROI_14 ROI_15 ROI_2 ROI_3 ROI_4
## 952 945 1315 1045 1190 1174 1497 1316 1155 1422 1097
## ROI_5 ROI_6 ROI_7 ROI_8 ROI_9
## 1420 1059 958 1114 945
##
## Cell types by ROI:
cross_tab <- table(codex_toy_data$ROI_num, codex_toy_data$Annotation5)
print(cross_tab[, c("BCL6- B Cell", "CD4 T", "CD8 T", "DC", "M1")])##
## BCL6- B Cell CD4 T CD8 T DC M1
## ROI_0 190 209 171 114 76
## ROI_1 189 208 170 113 76
## ROI_10 263 289 237 158 105
## ROI_11 209 230 188 125 84
## ROI_12 238 262 214 143 95
## ROI_13 235 258 211 141 94
## ROI_14 299 329 269 180 120
## ROI_15 263 290 237 158 105
## ROI_2 231 254 208 139 92
## ROI_3 284 313 256 171 114
## ROI_4 219 241 197 132 88
## ROI_5 284 312 256 170 114
## ROI_6 212 233 190 127 85
## ROI_7 191 211 172 115 77
## ROI_8 223 245 200 134 89
## ROI_9 189 208 170 113 76
BioGSP_Tutorial.R
# Color mapping for visualization
fill_colors <- c("#e07329", "#16964a", "#2958a8", "#dd2246", "#703594",
"#d0bd2a", "#69bd48", "#2297b1", "#b9e4f3", "#ae5c71", "#ab2170")
color_mapping <- data.frame(
Annotation5 = c("BCL6- B Cell", "BCL6+ B Cell", "CD4 T", "CD4 Treg", "CD8 T", "DC", "Endothelial", "M1", "M2", "Myeloid", "Other"),
Color = fill_colors
)
color_mapping_vector <- setNames(color_mapping$Color, color_mapping$Annotation5)
# Plot overall spatial distribution
p_overview <- ggplot(codex_toy_data, aes(x = X_cent, y = Y_cent, color = Annotation5)) +
geom_point(size = 0.8, alpha = 0.7) +
scale_color_manual(values = color_mapping_vector) +
labs(title = "Toy CODEX Spatial Cell Distribution",
x = "X Coordinate", y = "Y Coordinate") +
theme_minimal() +
scale_y_reverse() +
guides(color = guide_legend(override.aes = list(size = 2)))
print(p_overview)# Show distribution by ROI
p_by_roi <- ggplot(codex_toy_data, aes(x = X_cent, y = Y_cent, color = Annotation5)) +
geom_point(size = 0.8, alpha = 0.7) +
facet_wrap(~ROI_num, scales = "free") +
scale_color_manual(values = color_mapping_vector) +
labs(title = "Spatial Distribution by ROI",
x = "X Coordinate", y = "Y Coordinate") +
theme_minimal() +
scale_y_reverse() +
theme(legend.position = "bottom")
print(p_by_roi)# Highlight target cell types
target_data <- codex_toy_data %>% filter(Annotation5 %in% target_cells)
p_target <- ggplot(codex_toy_data, aes(x = X_cent, y = Y_cent)) +
geom_point(color = "lightgray", size = 0.3, alpha = 0.3) +
geom_point(data = target_data, aes(color = Annotation5), size = 1.2) +
facet_wrap(~ROI_num, scales = "free") +
scale_color_manual(values = color_mapping_vector[target_cells]) +
labs(title = "BCL6- B Cell and CD4 T Cell Distribution by ROI",
x = "X Coordinate", y = "Y Coordinate") +
theme_minimal() +
scale_y_reverse()
print(p_target)BioGSP_Tutorial.R
## Demonstrating different kernel family options available in SGWT...
# Compare different kernel families
kernel_comparison <- compare_kernel_families(x_range = c(0, 3),
scale_param = 1,
plot_results = TRUE)# Show available kernel family types
kernel_types <- c("mexican_hat", "meyer", "heat")
cat("\nAvailable kernel family types:\n")##
## Available kernel family types:
## 1. mexican_hat
## 2. meyer
## 3. heat
##
## Each kernel family provides:
## - mexican_hat: Gaussian scaling + LoG-style wavelet
## - meyer: Smooth cosine low-pass + band-pass transitions
## - heat: Exponential decay scaling + derivative-like wavelet
##
## Kernel values at different eigenvalues:
## x mexican_hat_scaling mexican_hat_wavelet meyer_scaling
## 1 0.00000000 1.0000000 0.0000000000 1
## 2 0.01507538 0.9998864 0.0002272412 1
## 3 0.03015075 0.9995456 0.0009086548 1
## 4 0.04522613 0.9989778 0.0020433121 1
## 5 0.06030151 0.9981835 0.0036296666 1
## 6 0.07537688 0.9971632 0.0056655569 1
## 7 0.09045226 0.9959176 0.0081482106 1
## 8 0.10552764 0.9944474 0.0110742486 1
## meyer_wavelet heat_scaling heat_wavelet
## 1 0 1.0000000 0.00000000
## 2 0 0.9850377 0.01484981
## 3 0 0.9702992 0.02925525
## 4 0 0.9557813 0.04322629
## 5 0 0.9414806 0.05677270
## 6 0 0.9273939 0.06990406
## 7 0 0.9135179 0.08262976
## 8 0 0.8998496 0.09495900
BioGSP_Tutorial.R
# Function to create binned data for SGWT analysis
create_binned_data <- function(cell_data, cell_type, x_bins = 30, y_bins = 30) {
# Filter for specific cell type
filtered_data <- cell_data %>% filter(Annotation5 == cell_type)
if (nrow(filtered_data) == 0) {
stop(paste("No cells found for cell type:", cell_type))
}
# Create bins
filtered_data <- filtered_data %>%
mutate(
x_bin = cut(X_cent, breaks = x_bins, labels = FALSE),
y_bin = cut(Y_cent, breaks = y_bins, labels = FALSE)
)
# Count cells in each bin
bin_counts <- filtered_data %>%
group_by(x_bin, y_bin) %>%
summarise(cell_count = n(), .groups = 'drop')
# Create complete grid
complete_grid <- expand.grid(x_bin = 1:x_bins, y_bin = 1:y_bins)
# Merge with counts (fill missing with 0)
binned_data <- complete_grid %>%
left_join(bin_counts, by = c("x_bin", "y_bin")) %>%
mutate(cell_count = ifelse(is.na(cell_count), 0, cell_count))
# Calculate proportions
binned_data$cell_proportion <- binned_data$cell_count / max(binned_data$cell_count)
return(binned_data)
}
# Create binned data for both cell types in ROI_1
# Filter for ROI_1 for focused analysis
roi1_data <- codex_toy_data %>% filter(ROI_num == "ROI_1")
bcl6nb_binned <- create_binned_data(roi1_data, "BCL6- B Cell")
cd4t_binned <- create_binned_data(roi1_data, "CD4 T")
cat("Binned data created successfully!\n")## Binned data created successfully!
## BCL6- B Cell bins with cells: 157 out of 900
## CD4 T bins with cells: 171 out of 900
# Visualize binned data
plot_binned_data <- function(binned_data, cell_type, color_val) {
ggplot(binned_data, aes(x = x_bin, y = y_bin, fill = cell_proportion)) +
geom_tile() +
scale_fill_gradient2(low = "white", high = color_val, name = "Proportion") +
labs(title = paste(cell_type, "Cell Distribution (30x30 bins)"),
x = "X Bin", y = "Y Bin") +
theme_minimal() +
scale_y_reverse() +
coord_equal()
}
p_bcl6nb_bins <- plot_binned_data(bcl6nb_binned, "BCL6- B Cell", color_mapping_vector["BCL6- B Cell"])
p_cd4t_bins <- plot_binned_data(cd4t_binned, "CD4 T", color_mapping_vector["CD4 T"])
# Show side by side
bins_combined <- p_bcl6nb_bins | p_cd4t_bins
print(bins_combined)BioGSP_Tutorial.R
# Prepare data for SGWT analysis
prepare_sgwt_data <- function(binned_data) {
sgwt_data <- binned_data %>%
mutate(
x = x_bin,
y = y_bin,
signal = cell_proportion
) %>%
select(x, y, signal)
return(sgwt_data)
}
# Prepare data for both cell types
bcl6nb_sgwt_data <- prepare_sgwt_data(bcl6nb_binned)
cd4t_sgwt_data <- prepare_sgwt_data(cd4t_binned)
cat("SGWT data prepared:\n")## SGWT data prepared:
## BCL6nB data points: 900
## CD4T data points: 900
## Signal range BCL6nB: 0 1
## Signal range CD4T: 0 1
BioGSP_Tutorial.R
## Applying SGWT to BCL6nB cell distribution...
# Apply SGWT to BCL6nB data
sgwt_bcl6nb <- SGWT(data.in = bcl6nb_sgwt_data,
signal = "signal",
k = 15, # Number of nearest neighbors
J = 3, # Number of scales
scaling_factor = 2,
kernel_type = "mexican_hat", # Unified kernel family
laplacian_type = "normalized",
k_fold = 2,
return_all = TRUE)## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.4886, 0.2443, 0.1221
## Performing SGWT decomposition...
## Reconstruction RMSE: 0.11396
## SGWT analysis completed for BCL6nB
## Reconstruction RMSE: 0.11396
##
## Applying SGWT to CD4T cell distribution...
sgwt_cd4t <- SGWT(data.in = cd4t_sgwt_data,
signal = "signal",
k = 15,
J = 3,
scaling_factor = 2,
kernel_type = "mexican_hat",
laplacian_type = "normalized",
k_fold = 2,
return_all = TRUE)## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.4886, 0.2443, 0.1221
## Performing SGWT decomposition...
## Reconstruction RMSE: 0.149172
## SGWT analysis completed for CD4T
## Reconstruction RMSE: 0.149172
# Energy analysis
energy_bcl6nb <- sgwt_energy_analysis(sgwt_bcl6nb)
energy_cd4t <- sgwt_energy_analysis(sgwt_cd4t)
# Combine energy data for comparison
energy_bcl6nb$cell_type <- "BCL6nB"
energy_cd4t$cell_type <- "CD4T"
energy_combined <- rbind(energy_bcl6nb, energy_cd4t)
# Plot energy distribution
p_energy <- ggplot(energy_combined, aes(x = scale, y = energy_ratio, fill = cell_type)) +
geom_bar(stat = "identity", position = "dodge") +
scale_fill_manual(values = c("BCL6nB" = color_mapping_vector["BCL6nB"],
"CD4T" = color_mapping_vector["CD4T"])) +
labs(title = "Energy Distribution Across SGWT Scales",
x = "Scale", y = "Energy Ratio", fill = "Cell Type") +
theme_minimal()
print(p_energy)BioGSP_Tutorial.R
## Testing different kernel families on BCL6nB data...
# Extract components for testing
eigenvalues <- sgwt_bcl6nb$graph_info$eigenvalues
eigenvectors <- sgwt_bcl6nb$graph_info$eigenvectors
signal_vector <- bcl6nb_sgwt_data$signal
# Test different kernel families
kernel_types_test <- c("mexican_hat", "meyer", "heat")
kernel_results <- list()
for (kernel_type in kernel_types_test) {
cat("Testing", kernel_type, "kernel family...\n")
result <- sgwt_forward(signal = signal_vector,
eigenvectors = eigenvectors,
eigenvalues = eigenvalues,
scales = sgwt_bcl6nb$parameters$scales,
kernel_type = kernel_type)
reconstructed <- sgwt_inverse(result)
rmse <- sqrt(mean((signal_vector - reconstructed)^2))
kernel_results[[kernel_type]] <- list(
result = result,
reconstructed = reconstructed,
rmse = rmse
)
cat("RMSE with", kernel_type, "kernel:", round(rmse, 6), "\n")
}## Testing mexican_hat kernel family...
## RMSE with mexican_hat kernel: 0.11396
## Testing meyer kernel family...
## RMSE with meyer kernel: 0.11527
## Testing heat kernel family...
## RMSE with heat kernel: 0.113988
# Create comparison table
kernel_comparison_table <- data.frame(
Kernel_Family = names(kernel_results),
RMSE = sapply(kernel_results, function(x) round(x$rmse, 6)),
Description = c("Gaussian scaling + LoG wavelet",
"Smooth cosine transitions",
"Exponential decay + derivative"),
stringsAsFactors = FALSE
)
kable(kernel_comparison_table,
caption = "Reconstruction Error Comparison Across Different Kernel Families")| Kernel_Family | RMSE | Description | |
|---|---|---|---|
| mexican_hat | mexican_hat | 0.113960 | Gaussian scaling + LoG wavelet |
| meyer | meyer | 0.115270 | Smooth cosine transitions |
| heat | heat | 0.113988 | Exponential decay + derivative |
BioGSP_Tutorial.R
# Function to visualize SGWT components
plot_sgwt_components <- function(sgwt_result, data_in, cell_type, color_val) {
coefficients <- sgwt_result$decomposition$coefficients
plot_data <- data_in
plots <- list()
# Original signal
plot_data$original <- sgwt_result$original_signal
p_orig <- ggplot(plot_data, aes(x = x, y = y, fill = original)) +
geom_tile() +
scale_fill_gradient2(low = "white", high = color_val) +
labs(title = paste(cell_type, "- Original")) +
theme_void() + scale_y_reverse() + coord_equal() +
theme(legend.position = "none")
plots[["original"]] <- p_orig
# Scaling function
plot_data$scaling <- Re(as.vector(coefficients[[1]]))
p_scaling <- ggplot(plot_data, aes(x = x, y = y, fill = scaling)) +
geom_tile() +
scale_fill_gradient2(low = "blue", mid = "white", high = "red") +
labs(title = "Scaling (Low-freq)") +
theme_void() + scale_y_reverse() + coord_equal() +
theme(legend.position = "none")
plots[["scaling"]] <- p_scaling
# Wavelet coefficients
for (i in 1:min(3, length(coefficients) - 1)) {
coeff_name <- paste0("wavelet_", i)
plot_data[[coeff_name]] <- Re(as.vector(coefficients[[i + 1]]))
p_wavelet <- ggplot(plot_data, aes_string(x = "x", y = "y", fill = coeff_name)) +
geom_tile() +
scale_fill_gradient2(low = "blue", mid = "white", high = "red") +
labs(title = paste("Wavelet", i)) +
theme_void() + scale_y_reverse() + coord_equal() +
theme(legend.position = "none")
plots[[coeff_name]] <- p_wavelet
}
# Reconstructed signal
plot_data$reconstructed <- sgwt_result$reconstructed_signal
p_recon <- ggplot(plot_data, aes(x = x, y = y, fill = reconstructed)) +
geom_tile() +
scale_fill_gradient2(low = "white", high = color_val) +
labs(title = "Reconstructed") +
theme_void() + scale_y_reverse() + coord_equal() +
theme(legend.position = "none")
plots[["reconstructed"]] <- p_recon
return(plots)
}
# Plot BCL6nB decomposition
bcl6nb_plots <- plot_sgwt_components(sgwt_bcl6nb, bcl6nb_sgwt_data, "BCL6nB",
color_mapping_vector["BCL6nB"])
bcl6nb_combined <- wrap_plots(bcl6nb_plots, ncol = 3)
print(bcl6nb_combined)# Plot CD4T decomposition
cd4t_plots <- plot_sgwt_components(sgwt_cd4t, cd4t_sgwt_data, "CD4T",
color_mapping_vector["CD4T"])
cd4t_combined <- wrap_plots(cd4t_plots, ncol = 3)
print(cd4t_combined)BioGSP_Tutorial.R
## Performing cross-correlation analysis between BCL6nB and CD4T...
# Calculate Graph Cross-Correlation (GCC)
gcc_result <- Cal_GCC(
data.in = data.frame(
BCL6nB = bcl6nb_sgwt_data$signal,
CD4T = cd4t_sgwt_data$signal
),
knee = c(length(eigenvalues) %/% 4),
signal1 = "BCL6nB",
signal2 = "CD4T",
eigenvector = eigenvectors
)
# Cosine similarity and Pearson correlation for comparison
cosine_sim <- cosine_similarity(bcl6nb_sgwt_data$signal, cd4t_sgwt_data$signal)
pearson_cor <- cor(bcl6nb_sgwt_data$signal, cd4t_sgwt_data$signal)
cat("Graph Cross-Correlation between BCL6nB and CD4T:", round(gcc_result, 4), "\n")## Graph Cross-Correlation between BCL6nB and CD4T: 0.9412
## Cosine similarity between BCL6nB and CD4T: 0.3414
## Pearson correlation between BCL6nB and CD4T: 0.2188
# Create comparison table
correlation_comparison <- data.frame(
Method = c("Graph Cross-Correlation (GCC)", "Cosine Similarity", "Pearson Correlation"),
Value = c(round(gcc_result, 4), round(cosine_sim, 4), round(pearson_cor, 4)),
Description = c("Frequency-domain correlation using graph structure",
"Geometric similarity measure",
"Linear correlation coefficient")
)
kable(correlation_comparison,
caption = "Cross-Cell Type Correlation Analysis: BCL6nB vs CD4T")| Method | Value | Description |
|---|---|---|
| Graph Cross-Correlation (GCC) | 0.9412 | Frequency-domain correlation using graph structure |
| Cosine Similarity | 0.3414 | Geometric similarity measure |
| Pearson Correlation | 0.2188 | Linear correlation coefficient |
BioGSP_Tutorial.R
## === Demonstrating Cross-ROI Analysis ===
# Compare BCL6- B Cell patterns across different ROIs
roi_comparison <- list()
for (roi in c("ROI_1", "ROI_2", "ROI_3")) {
roi_data <- codex_toy_data %>% filter(ROI_num == roi)
roi_bcl6nb <- roi_data %>% filter(Annotation5 == "BCL6- B Cell")
cat("ROI", roi, "- BCL6- B Cell:", nrow(roi_bcl6nb), "\n")
cat(" X range:", round(range(roi_bcl6nb$X_cent), 2), "\n")
cat(" Y range:", round(range(roi_bcl6nb$Y_cent), 2), "\n")
# Quick binning for comparison
if (nrow(roi_bcl6nb) > 20) { # Only if sufficient cells
binned <- create_binned_data(roi_data, "BCL6- B Cell", x_bins = 15, y_bins = 15)
sgwt_data <- prepare_sgwt_data(binned)
# Quick SGWT analysis
result <- SGWT(data.in = sgwt_data,
signal = "signal",
k = min(8, nrow(sgwt_data) - 1),
J = 2, # Fewer scales for quick analysis
kernel_type = "mexican_hat",
k_fold = 2,
return_all = TRUE)
roi_comparison[[roi]] <- list(
cells = nrow(roi_bcl6nb),
rmse = result$reconstruction_error,
active_bins = sum(binned$cell_count > 0)
)
cat(" SGWT RMSE:", round(result$reconstruction_error, 4), "\n")
cat(" Active bins:", sum(binned$cell_count > 0), "out of 225\n\n")
}
}## ROI ROI_1 - BCL6- B Cell: 189
## X range: 28.14 59.85
## Y range: 76.06 106.72
## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.5156, 0.2578
## Performing SGWT decomposition...
## Reconstruction RMSE: 0.1359
## SGWT RMSE: 0.1359
## Active bins: 109 out of 225
##
## ROI ROI_2 - BCL6- B Cell: 231
## X range: 49.92 87.46
## Y range: 13.05 46.97
## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.5156, 0.2578
## Performing SGWT decomposition...
## Reconstruction RMSE: 0.099712
## SGWT RMSE: 0.0997
## Active bins: 99 out of 225
##
## ROI ROI_3 - BCL6- B Cell: 284
## X range: 36.21 95.79
## Y range: 16.7 70.47
## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.5156, 0.2578
## Performing SGWT decomposition...
## Reconstruction RMSE: 0.074979
## SGWT RMSE: 0.075
## Active bins: 96 out of 225
# Create comparison table
if (length(roi_comparison) > 0) {
comparison_df <- data.frame(
ROI = names(roi_comparison),
BCL6_B_Cells = sapply(roi_comparison, function(x) x$cells),
Active_Bins = sapply(roi_comparison, function(x) x$active_bins),
SGWT_RMSE = round(sapply(roi_comparison, function(x) x$rmse), 4)
)
kable(comparison_df, caption = "Cross-ROI BCL6- B Cell Analysis Comparison")
}| ROI | BCL6_B_Cells | Active_Bins | SGWT_RMSE | |
|---|---|---|---|---|
| ROI_1 | ROI_1 | 189 | 109 | 0.1359 |
| ROI_2 | ROI_2 | 231 | 99 | 0.0997 |
| ROI_3 | ROI_3 | 284 | 96 | 0.0750 |
## Cross-ROI analysis demonstrates how spatial patterns vary across regions!
BioGSP_Tutorial.R
## === SGWT Analysis Summary ===
## Dataset: Multi-ROI Toy CODEX spatial data
## - Total cells analyzed: 18604
## - Number of ROIs: 16
## - ROI_1 focused analysis (30x30 grid): 900 bins
## - BCL6nB active bins (ROI_1): 157
## - CD4T active bins (ROI_1): 171
##
## SGWT Analysis Results:
## - BCL6nB reconstruction RMSE: 0.11396
## - CD4T reconstruction RMSE: 0.149172
## - Number of scales: 3
## - Kernel type: mexican_hat
##
## Cross-Correlation Results:
## - Graph Cross-Correlation: 0.9412
## - Cosine Similarity: 0.3414
## - Pearson Correlation: 0.2188
if (abs(gcc_result) > 0.3) {
cat("\n✓ Moderate to strong spatial correlation detected between cell types\n")
} else {
cat("\n✓ Weak spatial correlation between cell types\n")
}##
## ✓ Moderate to strong spatial correlation detected between cell types
##
## ✓ Tutorial completed successfully using the BioGSP package!
BioGSP_Tutorial.R
This tutorial demonstrated the key capabilities of the BioGSP package:
The BioGSP package provides a comprehensive framework for spatial analysis of multiplexed imaging data, enabling researchers to uncover multi-scale spatial patterns and relationships in complex tissue environments.
For more information, see ?codex_toy_data for dataset details and explore other BioGSP functions using help(package = "BioGSP").